function [V_r,A_polr,W_polr] = vret_4(A,W,last,a,b_db,b_ss,optfunc,solver_param,taxf,mortp_80)
   
    sigma  = solver_param(1);
    beta   = solver_param(2);
    aT     = solver_param(3);
    aS     = solver_param(4);
    r      = solver_param(5);
    V_term = solver_param(6);
    scaler = solver_param(7);
    pi     = solver_param(8);
    mortp  = mortp_80(aT-aS-a);
    itau    = taxf{3,1}; 
    tau_k  = 0;
    
	%Inputs:
	%A = state space for assets
	%W = state space for DC balance
	%f = state space for labor disutility
	%last = 1 if computing last period functions
	%itau = tax rate (or function)
    %tau_k = tax on capital     
	%r = interest rate
	%aT = Terminal age
	%aS = start age
	%a = current age
	%b_db = DB annuity payout 
    %b_ss = SS annuity payout
    %scaler = consumption utililty scaler
	%optfunc = cell that holds previously estimated value and policy functions
	%V_term = terminal value (=0)
	
	%Outputs 
	%V_r = value function when retired 
	%A_polr = policy function for saving non-DC assets
	%W_polr = policy function for saving DC assets
    
    %Set tax rates
    if last == 0
        tau = itau((aT-aS-a-1)*(length(W)^2)+1:(aT-aS-a)*length(W)^2);
        tau = reshape(tau,[length(W),length(W)]);
    elseif last == 1
        tau = itau((aT-aS-1)*(length(W)^2)+1:(aT-aS)*length(W)^2);
        tau = reshape(tau,[length(W),length(W)]);
    end   
    

    b_ss62 = b_ss(13); %claim age 62 SS annuity
    %b_db and b_ss are the vectors of potential annuity levels (initial age to current
    %age)
    b_db = b_db(1:aT-aS-a);
    b_ss = b_ss(1:aT-aS-a);
    
    %Cells to hold policy and value functions
    %V_r = value function when retired 
	%A_polr = policy function for saving non-DC assets
	%W_polr = policy function for saving DC assets
    V_r = cell(aT-aS,1);
    A_polr = cell(aT-aS,1);
    W_polr = cell(aT-aS,1);
    
    %Pre-create state space grid 
    W_0 = W.*(1+r)^(aT-aS-a-1); % DC balance grid current period
    W_1 = W.*(1+r)^(aT-aS-a);   % DC balance grid next period
    A_0 = A.*(1+r)^(aT-aS-a-1); % Non-DC grid current period
    A_1 = A.*(1+r)^(aT-aS-a);   % Non-DC grid next period
    [A_nd,W_nd] = ndgrid(A_1,W_1);

        for bb = 1:length(b_db) %loop over b space (1 is youngest, length(b) is oldest)
            %Pre-create empty value and policy functions for
            %each loop iteration
            V_rT = zeros(size(A_nd));
            A_polrT = zeros(size(A_nd));
            W_polrT = zeros(size(A_nd)); 
            for i = 1:length(A_0) %loop over A space 
                for j = 1:length(W_0) %loop over W space   

                        %Pick bt (annuity level: for DB this is deflated from ret age to current age)
                        if bb < 13 && a > 18      % retired before 62 and current age <62: get b_DB(R_age) 
                            bt = b_db(bb)*(1/(1+pi))^(aT-aS-bb-a);                     
                        elseif bb < 13 && a <= 18      % retired before 62 and current age >=62: get b_DB(R_age) + b_SS62
                            bt = b_db(bb)*(1/(1+pi))^(aT-aS-bb-a) +  b_ss62;
                        elseif bb > 13                       % retired after 62 and any current age: get b_DB(R_age) + b_SS(R_age)
                            bt = b_db(bb)*(1/(1+pi))^(aT-aS-bb-a) + b_ss(bb);
                        end
                        %Pick continuation value
                        if last == 1
                            V_cont = V_term;
                        elseif last == 0
                            V_cont1 = optfunc{a,4}; %bb'th element of the cell
                            V_cont  = V_cont1{bb};
                        end

                        %Determine consumption (state depdendent on whether
                        %drawing down DC wealth)

                        %Distributions from DC accounts
                        distr = W_0(j)-W_nd./(1+r);     %distribution from DC assets
                        T=(W_0(j)-W_nd./(1+r)>=0);      %negative distributions
                        distr = distr.*T;             %zero out negative distributions
                        %Deduct taxes
                        totrinc = bt+distr; %total retirement income from pensions + SS
                        for k = 1:length(A_0)
                            totrinc(k,:) = totrinc(k,:).*(1-tau(j,:));
                        end
                        %Consumption net of taxes
                        cons = totrinc+A_0(i)-A_nd./(1+r*(1-tau_k));
                        cons(cons<0)=0; %zero out negative consumption levels

                        %Value function
                        v_r_temp = scaler.*(1/(1-(1/sigma))).*(cons).^(1-(1/sigma))...
                           +(1-mortp).*beta.*V_cont+500;
                        %Disqualify building up W when retired
                        v_r_temp(:,j+1:end) = -Inf;
                        %Eliminate neg consumption levels from maxima
                        T = (cons>0);
                        v_r_temp = T.*v_r_temp + (1-T).*-1e50;

                        %Determine maxima 
                        %max over W dimension
                        [dimwmax,ind1]=max(v_r_temp,[],2);
                        %max over A dimension
                        [vmax,ind2]=max(dimwmax);
                        %index of W policy (pick ind2th row from ind1th column)
                        Wind = ind1(ind2);       
                        %index of A policy
                        Aind = ind2;

                        %Assign value and policy functions
                        V_rT(i,j)=vmax; 
                        A_polrT(i,j)=A_1(Aind);
                        W_polrT(i,j)=W_1(Wind);
                end
            end
            
            V_r{bb} = V_rT;
            A_polr{bb} = A_polrT;
            W_polr{bb} = W_polrT;
            
        end
	
   
   end
